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ABSTRACT 

Matter mixing is one important topic in the study of core-collapse supernova 
(CCSN) explosions. In this paper, we perform two-dimensional hydrodynamic 
simulations to reproduce the high velocity 56 Ni clumps observed in SN 1987A. 

This is the first time that large density perturbation is proposed in the CCSN 
progenitor to generate Rayleigh-Taylor (RT) instability and make the effective 
matter mixing. In the case of a spherical explosion, RT instability is efficient 
at both C+O/He and He/H interfaces of the SN progenitor. Radial coherent 
structures shown in perturbation patterns are important for obtaining high ve¬ 
locity 56 Ni clumps. We can also obtain matter mixing features and high velocity 
56 Ni clumps in some cases of aspherical explosion. We find that one of the most 
favorable models in our work has a combination of bipolar and equatorially asym¬ 
metric explosions in which at least 25% of density perturbation is introduced at 
different composition interfaces of the CCSN progenitor. These simulation results 
are comparable to the observational findings of SN 1987A. 

Subject headings: Hydrodynamics — instabilities — nuclear reactions, nucleosyn¬ 
thesis, abundances — shock waves — supernovae: general 
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Introduction 


Heavy elements in the universe can be released by core-collapse supernova (CCSN) 
explosions. It has been found in the observations of SN 1987A that some fast-moving 56 Ni 
clumps ejected outside have similar velocities of helium and hydrogen (e.g., Hanuschik et 
al. 1988). This phenomenon is called matter mixing in this paper. The radioactive decay 
of 56 Ni—)• 56 Co—)• 56 Fe (the half-lives of 56 Ni and 56 Co are 5.9 and 77.1 days, respectively) 
and the related energy heating of SN 1987A provide excellent opportunities to study matter 


mixing problems and explosion morphology. Early detections of hard X-ray emission (Dotani 


et al. 


1987| Sunyaev et al. 1987) and gamma-ray lines from 56 Co (Matz et al. 1988 Leising 


& Share 1990) in SN 1987A indicate a high velocity component of radioactive 56 Ni from the 


CCSN explosive nucleosynthesis. The explosion energy of SN 1987A is given as E/M env = 
(1.1 — 1.5) x 10 50 erg M 0 4 where M env is the mass of the hydrogen-rich envelope in the 
supernova (SN) progenitor and the total mass of synthesized 56 Ni in SN 1987A is about 
O.O7M 0 (Shigeyama et al. 1988; Shigeyama & Nornoto 1990; see also Handy et al. 2014 for 
the summary of SN 1987A explosion energy). The mass of high velocity 56 Ni clumps was 


determined to be about a few times of 10 3 M 0 ( 

Utrobin et al.||l995 

) from the measurement 

of Ha fine structure ( 

Hanuschik et al. 

1988). Haas et al. 

(1990 

) and 

Spyromilio et al. 

(1990) 


observed doppler-broadened line profiles of [Fe II] in SN 1987A. Haas et al. (1990) found 
that the velocity can reach about 4,000 km s -1 , and 4 — 17% of the total 56 Ni has a high 
velocity larger than 3,000 km s -1 . Furthermore, the line profile is asymmetric, and its peak is 
located in the redshifted side of about 1,000 km s -1 (Haas et al.|l990| Spyromilio et al. 1990). 


Another interesting radioactive element produced by the explosive nucleosynthesis is 44 Ti. 
Models of the spectrum and lightcurve of SN 1987A powered by 44 Ti have been presented 
( Fransson fe Kozma|l993 Jerkstrand et al.| 2011). Although the half-life of 44 Ti—)- 44 Sc—)- 44 Ca 


is about 59 years (the half-lives of 44 Ti and 44 Sc are 58.9 ±0.3 years and 3.97 hr, respectively, 
see Ahmad et al. 2006), hard X-ray emission lines from the decay of 44 Ti in SN 1987A were 
also detected and the estimated mass of 44 Ti is (3.1 ±0.8) x 10 -4 M 0 
The inner ejecta of SN 1987A is heated by the radioactivity of 44 


(Grebenev et al. 

2012) 

Ti ( 

Kjaer et al. 

2010) 


Moreover, three-dimensional (3D) structures reconstructed by observations show that the 


inner ejecta of SN 1987A may lie in the plane of the equatorial inner ring (Kjaer et al. 2010 


Larsson et al. 2013 Sinnott et al. 2013). As the collision between this inner ejecta and the 


edge of inner ring may happen, this morphology can naturally explain the observed inner 
ejecta rebrightening (e.g., Groningsson et al. 2008, Dwek et al. 2010). 

The physical mechanisms of CCSN explosions are under debate. In general, there are at 
least two interesting aspects of CCSN explosion mechanisms (e.g., Janka 2012; Kotake et al. 
2012; Burrows 2013). In the framework of a neutrino heating mechanism, standing accretion 


shock instability (SASI, Blondin et al. 2003) and neutrino-driven convection (Herant et al. 
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1994 Burrows et al. 1995) were proposed to generate a globally anisotropic structure and 


enhance the effective heating. Many hydrodynamic simulations have been performed during 
recent years (e.g., Nordhaus et al. 2010; Suwa et al. 2011; Kuroda et al. 2012; Miiller et 
al. 2012; Bruenn et al. 2013; Couch & Ott 2013; Hanke et al. 2013; Murphy et al. 2013; 
Ott et al. 2013, Wongwathanarat et al. 2013; Takiwaki et al. 2014). On the other hand, 
MHD effects can trigger a jet-like CCSN explosion, as shown in some numerical simulations 


(Kotake et al. 2004; Sawai et al. 2005, 2013 Burrows et al. 2007; Takiwaki et al. 2009). 


So far, some tests to solve the matter mixing problems of CCSN explosions have been 
carried out under the hydrodynamic framework, and the explosive nucleosynthesis has been 
adopted in particular. Raylcigh-Taylor (RT) instability (or Richtmyer-Meshkov instability 
in the condition of entropy impulsive acceleration) can happen naturally in CCSN explosion 
models, and it plays a vital role on shock propagation and matter mixing. For example, 
it was found in early works that He/H and C+O/He layers in an SN progenitor can be 
unstable against strong RT instability ( Ebisuzaki et al.||1989 Benz fe Thielemann||1990 ). RT 
fingers first dominate at the H/He interface, and at late times are prominently shown at the 


He/C+O interface (Muller et al. 1991). 


Perturbation in the progenitor star has been considered to be an origin of RT insta¬ 
bility and this mechanism may be adopted to solve matter mixing problems. In principle, 
perturbations are related to the physics of convection and turbulence in the stellar interior. 
The integrated convection turnover time across the whole convective region is comparable 
to the nuclear timescale. This indicates that the composition in the CCSN progenitor is not 


homogeneous, thus the progenitor system has strong fluctuations (Bazan & Arnett 1998). 


Stellar convection is thought to make significant perturbations based on the dynamical sim¬ 
ulations of progenitor stars (e.g., Bazan & Arnett 1994; Bazan & Arnett 1998; Meakin & 


Arnett 2007b; Arnett & Meakin 2011a). Bazan & Arnett (1994) suggested progenitor per 


turbations of density, temperature, pressure, and electron fraction in the order of 5%, 1%, 


3%, and 0.08%, respectively. Arnett (1994) expected 2.5% of density perturbation and 0. 
of velocity perturbation in the oxygen shell of progenitor star. CCSN progenitor features 
were also explored by an extended work of density perturbation at the edges of metal con¬ 


vective zones (Bazan & Arnett 1998). Furthermore, some violent dynamic eruptions in the 


convective zones with simultaneous carbon and oxygen shell burning inside the progenitor 
star were found in 3D simulations, and about 10% of density perturbation was detected in 


those convective zones (Meakin & Arnett 2007b; Arnett & Meakin 2011a). 


The position that RT instability occurs in the CCSN progenitor has also been discussed. 


Arnett et al. (1989) noticed the possibility that RT instability begins at the edge of the dense 


mass shell, where heavy elements mix with hydrogen through the He envelope. Muller et al. 
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(1991) mentioned that RT fingers form at both He/C+O and H/He interfaces. Nagataki et 


al. (1998) introduced velocity perturbation at the layer between ffe and ff. Ellinger et ah 


(2012) found that this instability at He/H interface is caused by reverse shock. Hammer et 


al. (2010) and Kifonidis et al. (2003, 2006) found an instability at the Si/O interface, which 


was not shown in the results of Ellinger et al. (2012). Although all of instabilities at Si/O, 


He/C+O, and He/H interfaces were described by Kifonidis et al. (2003), the instability at the 
C+O/He boundary may be the most prominent one (Ellinger et al. 2012). This is consistent 
with the conclusion given by Kane et al. (2000). 


It is expected that perturbations in the CCSN progenitor may generate or enhance the 
matter mixing. Simulations with introducing an initial density perturbation of 5% amplitude 
in the progenitor model failed to reproduce the high velocity (> 3500 km s -1 ) 56 Ni clumps 
(Herant & Benz 1991, 1992). On the other hand, the possibility of large perturbations 


in the CCSN progenitor has been investigated. The matter mixing can be reproduced by 
simulation models if 30% of the initial velocity perturbation and 30% of the initial energy 


perturbation are introduced in the CCSN progenitor (Yamada & Sato 1991). Nagataki et 


al. (1998) introduced 30% of velocity perturbation at the H/He boundary layer in the CCSN 
progenitor and obtained high velocity (> 3000 km s -1 ) 56 Ni clumps. However, the spacial 
resolution of these simulations is rather low. Three-dimensional simulations with different 


metallicities in the CCSN progenitor were also performed to develop RT instability (Joggerst 
2010a|b). It is interesting that large asymmetric structures during late times 


et al. 2009 


of a CCSN explosion may be originated by the significant perturbation that arises from an 
initial explosion itself, but not from merging RT structures. 

The asymmetric effect of a CCSN explosion is also involved as an important point for the 


matter mixing. Nagataki et al. (1998) and Nagataki (2000) found that asymmetric or jet-like 


explosions can effectively change the velocity distribution of 56 Ni. Hungerford et al. (2003 


2005) performed 3D smoothed particle hydrodynamics (SPH) simulations and compared the 


difference between symmetric and asymmetric CCSN explosions. They confirmed that an 
asymmetric explosion makes the matter mixing in supernova ejecta, although high velocity 


56 Ni clumps could not be reproduced. Kifonidis et al. (2003) comprehensively discussed RT 


instability and explosion evolution in one-dimensional (ID) and two-dimensional (2D) cases. 
Introducing low-mode perturbation in a CCSN explosion progenitor may reproduce high 
velocity (> 3300 km s -1 ) 56 Ni clumps (Kifonidis et al. 2006), although robust conclusions 


of 2D simulations were further discussed by Gawryszczak et al. (2010). Three-dimensional 


simulations of CCSN explosions were also performed by Hammer et al. (2010), and high 
velocity 56 Ni could be reproduced. However, their models neglected the effect of matter 
fallback into the compact remnant. As recently suggested by Ono et al. (2013), not only 
perturbations in the SN progenitor but also asymmetric explosions (clumpy morphologies 
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in explosions) may be necessary to reproduce the observed high velocity 56 Ni clumps of SN 
1987A. 

In this paper we focus on RT instability as one possible solution of matter mixing prob¬ 
lems in SN 1987A. We have already performed 2D hydrodynamic simulations to examine 
the matter mixing of CCSN (see Paper I). If a large radial velocity perturbation of 30% 
amplitude is introduced just before the shock wave reaches the hydrogen envelope of the 
progenitor star, high velocity 56 Ni clumps of 3,000 km s -1 may be reproduced, in the case 
of bipolar explosion. This indicates that RT instability originating from large perturbations 
in the progenitor may be one reason for the matter mixing in a CCSN explosion. Here, we 
further extend the work in Paper I. In particular, this is the first time that large density 
perturbation is systematically considered in a CCSN progenitor to make strong RT instabil¬ 
ity, and consequently dig out a part of inner heavy elements to outer layers during CCSN 
explosion. Here, we mention some possible motivations of the large density perturbations 
introduced in our CCSN progenitor. The global turbulent feature of the CCSN progenitor 
interior was explored ( Meakin fe Arnett|2007a Arnett et al. [2009 Arnett fe Meakin||2011a|b 


Viallet et al. 2013), and it was realized that a realistic progenitor model is crucial for CCSN 
explosion (Smith & ArnettS2014). In this paper, we adopt the CCSN progenitor that has been 
successfully used to produce the lightcurve, 56 Ni mass, and explosion energy of SN 1987A 


(Shigeyama et al. 1987, 1988 Hashimoto et al. 1989 Shigeyama & Nomoto 1990). Then, 


we introduce large density perturbations to this progenitor. We suggest that this progenitor 
with a large density perturbation may be physically equivalent/corresponding to a realistic 
model with a global turbulent feature and strong dynamic instabilities. Furthermore, Figure 
2 shows that the radial coherent features at large scales exist in some of our perturbation 
patterns. These radial coherent features are important for reproducing the matter mixing of 
SN 1987A. The coherent features shown in our perturbation patterns may capture the similar 


behavior driven by the turbulent convection in the stellar interior (Meakin & Arnett 2007a 


Arnett et al. 2009 

Arnett & Meakin 2011a|b 

Viallet et al. 2013 

). Moreover, 

Smith & Arnett 


(2014) argued two important effects of the hydrodynamic instability and turbulence in the 


progenitor model: (1) finite amplitude fluctuations in velocity and temperature, and (2) 
nonlinear interaction with nuclear burning during the latest stages of the stellar evolution. 


These effects are not included in current stellar evolution models. Smith (2007) claimed that 


SN 1987A is related to the luminous blue variable (LBV), which shows violent dynamical 
eruption or enhanced mass loss prior to the core-collapse. This invokes large density pertur¬ 
bation in the hydrogen envelope of the CCSN progenitor. We also propose Betalgeuse as an 
observational example: The clumpy structure in the inner envelope may be related to the 
giant convection cells of the outer atmosphere (Decin et al. 2012). The convection pattern 


can be also unveiled (Montarges et al. 2014). These observational evidences also suggest 
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large density perturbation in the hydrogen envelope of the CCSN progenitor. 

We suggest some physical factors—such as C+O layer convection, non-radial pulsa¬ 
tions in LBV stars, cool hydrogen envelopes, and convection-driven asymmetries in the deep 
interior— are affecting the development of RT instability and the matter mixing of SN 1987A. 
Although we introduce large density perturbations in the progenitor model as a major cause 
to develop RT instability and produce the matter mixing in this paper, there are many path¬ 
ways to the final asymmetries after the shock breakout. We focus on perturbations in the 
progenitor star, and its growth clue to RT instability during the shock wave propagation. 
As mentioned, perturbations could be introduced by combinations of different sources, such 
as convection in very deep, cool hydrogen envelopes; convective shell burnings in the deep 
interior; and nonlinear, non-radial pulsations in an LBV star. Aside from the origins of the 
perturbations, the growth of RT instability is affected by the progenitor structure. For giant 
RT fingers, a large hydrogen envelope, as in normal red supergiant (RSG) stars, is required 
(e.g., Wongwathanarat et al. 2015). On the other hand, RT instabilities in Wolf-Rayet pro¬ 
genitors may be significantly reduced by the lack of an H/He boundary with a strong density 
gradient. However, deeper asymmetries stemming from the RT instabilities at inner com¬ 
position boundaries and/or an aspherical explosion could still produce structure. Because 
the progenitor of SN 1987A is known to be a blue supergiant (BSG), which had a compact 
hydrogen envelope, it is difficult to produce prominent RT fingers compared to those of RSG 
stars. In addition, to dig out the innermost ejecta to outer layers, other factors-such as the 
compactness of the C+O shell, the density structure of He shell, and the density gradient at 
He/H boundary-may also be important (Wongwathanarat et al. 2015). 

We propose a large density perturbation in the progenitor model as follows: We test 
three cases of perturbations for spherical explosion. For all cases, the amplitude of pertur¬ 
bations are randomly assigned, but the length scale of the perturbations are different. First, 
perturbations are introduced at the interfaces of different composition layers, with an expo¬ 
nential decay of its amplitude in the radial direction. At each boundary layer, the length 
scale of the perturbation in the angular direction is controlled by a fixed wave number. In 
this case, perturbations have the radial coherent structure. Second, the perturbation length 
scale in the radial direction is proportional to the density scale height determined by the 
progenitor model. The perturbation length scale in angular direction is controlled by a fixed 
wave number. Third, we constrain the perturbation length scale in both radial direction 
and angular direction, as it is proportional to the density scale hight. All of these pertur¬ 
bation modes are comprehensively described in Section 2.3 and illustrated in Figures 1 and 
2. We also test three cases for aspherical explosion under the condition that perturbations 
are introduced at interfaces of different composition layers. First, we take the case that 
explosion along the polar direction is stronger than that along the equatorial direction. This 
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corresponds to a bipolar or jet-like SN explosion. Second, we propose the case of an SN 
explosion that is asymmetric to the equatorial plane. Third, we combine the two cases to 
take into account both bipolar and equatorially asymmetric effects of the CCSN explosion. 
The observed [Fe 11] line in SN 1987A shows an asymmetric feature in the velocity profile 


(Haas et al. 1990 Spyromilio et al. 

1990 

). The observation of SN 1987A in 3D structure 

presents an elongated morphology ( 

\jaer et al. 

2010 

). These observational evidences invoke 


us to make a combination of mildly bipolar and equatorially asymmetric explosions. In our 
work, we concern the radial velocity profile and line of sight velocity profile of 56 Ni, the total 
ejected mass of 56 Ni, and the mass of high velocity (> 3000 km s -1 ) 56 Ni. We also discuss 
the total ejected mass of 44 Ti and the related explosion energy. In Section 2, we review the 
procedures of the numerical simulation presented in Paper I. Numerical simulation methods 
and the progenitor star are briefly described and many density perturbation modes are given 
in detail. In Section 3, we comprehensively show the matter mixing results of 56 Ni in our 
simulations. Some important issues associated with our results are discussed in Section 4. 
Finally, we summarize our achievements in Section 5. 


2. Numerical Simulations 

We now briefly describe numerical methods, initial conditions, and progenitor model, 
which are presented in Paper I. We list several key points for our simulation models in 
Section 2.1 and 2.2. We extend the work of Paper I to explore the possibility of reproduc¬ 
ing high velocity 56 Ni clumps by using the progenitor with large density perturbation, and 
comprehensive perturbation modes are presented in Section 2.3. 


2.1. Numerical Arrangements and Initial Conditions 


We perform an adaptive mesh refinement (AMR) hydrodynamic code, FLASH (Fryxell 
et al. 2000), and a 2D spherical coordinate (r, 6) is adopted. The initial computation domain 
covers the region with a radius of 1.4 x 10 s cm < r < 3.0 x 10 9 cm and an angle of 0 < 9 < tt. 
The outer boundary of the progenitor model is limited at r = 3.4 x 10 12 cm, which corresponds 
to the surface of the hydrogen-rich envelope. A wind component is added beyond this radius. 
The base-level grid points are 48(r) x 12(0) and the maximum refinement level is set to be 
7. Thus, the effective maximum numbers of grid points are 3072(r) x 768(0). The minimum 
effective cell sizes are about 10 km and 0.23 deg in radial direction and angular direction, 
respectively. As the SN shock propagates and reaches close to the radial outer boundary, the 
computation domain is extended via remappings. When the physical values of the extended 














regions are set to the values of the progenitor model, 44 remappings are required. We use 
the monotonic cubic interpolation scheme to interpolate physical values (Steffen 1990). The 
computational cost is about 20 thousand CPU hr for each simulation model. 


The remapping and boundary conditions applied in this work were presented in Paper 
I. If the forward shock and/or RT fingers reach close to the radial outer boundary (approxi¬ 
mately 20% of the radial computational domain from the radial outer boundary), the radial 
size of the computational region is extended by a factor of 1.2. At the radial inner boundary, 
a “reflection 1 ' boundary condition is employed at the start of simulation. After the forward 
shock has reached the composition interface of C+O/He (corresponding to the radius of 
6 x 10 9 cm), it is switched to a “diode” boundary condition that allows matter to flow out of 
the computational domain, but inhibits it from entering the computational domain through 
the inner boundary in order to include the effect of matter fallback. Although changing the 
switch timing can somewhat affect the degree of the innermost matter fallback, we fix the 
switch timing to make sure that the mass of 56 Ni that remains in the computational domain 
does not become too small compared with that for SN 1987A. At the radial outer boundary, 
we adopt the ’’diode” boundary condition over the simulation. In the angular direction, we 
take the “reflection” boundary condition. Meanwhile, the mass inside the inner computation 
domain is regarded as a point source, and both the point source gravity and spherically aver¬ 
aged self-gravity are considered. The Helmholtz equation of state (EOS) is adopted for our 
simulations, whereas a modified ideal gas and radiation EOS is used for a very low density 
region, and we blend the two kinds of EOS at the transition density region (see Paper 1 for 
the details). 

Explosive nucleosynthesis is applied and we obtain the following elements from FLASH 
code: n, p, 1 H, 3 He, 4 He, 12 C, 14 N, 16 0, 20 Ne, 24 Mg, 28 Si, 32 S, 36 Ar, 40 Ca, 44 Ti, 48 Cr, 52 Fe, 
54 Fe, and 56 Ni. Advection equations are solved to trace the distribution of these elements. 
Energy depositions due to the radioactive decay of 56 Ni to 56 Fe are included, as described in 
Paper 1. However, some isotopes of iron group elements (e.g., 5 ' Ni and 58 Ni) are not included 
in the original FLASH code. Therefore, the masses of 56 Ni and 44 Ti obtained by FLASH code 
may be overestimated. In order to clarify the effect of this limited nuclear reaction network 
on the mass of heavy elements, an explosive nucleosynthesis calculation with a large nuclear 
network is performed as a post-process. We stress this issue in Section 4. The process is 
presented in the Appendix. 
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2.2. Progenitor Model 


The progenitor model of CCSN explosion adopted in this paper is a 16.3 M 0 BSG 
star with 6 M 0 helium core and 10.3 M 0 hydrogen envelope. This progenitor model has a 
radius set as 3.4 x 10 12 cm, and has been successfully adopted for reproducing the lightcurve, 
explosion energy, and total 56 Ni mass (0.07 M 0 ) of SN 1987A (Shigeyama et al. 1987, 1988 


Nomoto & Hashimoto 

1988 

Shigeyama & Nomoto 

1990) 


A wind component with the 
density profile of p oc r~ 2 is added and the uniform temperature of T = 10 4 K is adopted. 
The inner density of the wind component is 3.0 x 10~ 10 g cm -3 . The wind component is 
extended to the radius of 4.5 x 10 12 cm and simulations are carried out until just before the 
shock waves reach the radius. The density of the wind component is smoothly connected to 
that of the stellar surface. If we assume that the wind velocity is 15 km s _1 , the mass loss 
rate is 1.8 x 10 -3 M 0 yr _1 ; this value is larger than 1.0 x 10~ 5 M 0 yr^ 1 used in Gawryszczak 
et al. (2010). We restrict our simulation to just after the shock breakout, and the selected 
wind profile does not affect the matter mixing. 


We set the total input energy to be 2.5 x 10 51 erg in most simulation models and 
divide this input energy into two equal parts: kinetic energy and thermal energy. The final 
explosion energy is related to the input energy. Modification of total input energy value is 
also attempted in the two simulation models. 


Si, C+O, He, and H layers are settled as an onion-like structure in this progenitor. The 
outer boundary of the Si-rich layer is at the radius of 3 x 10 s cm, the interface between C+O 
and He core is at 6 x 10 9 cm, and the hydrogen envelope begins at 5 x 10 10 cm. It is also 
noted that oxygen shell burning is limited at 3 x 10 9 cm 

Therefore, we propose that these numbers are four critical radii in the progenitor star. 


(Arnett|l994 Bazan & Arnett 1998). 


2.3. Explosion and Density Perturbation Modes 

2.3.1. Spherical Explosion of SCI, SC2, and SC3 


In this subsection, we introduce three density perturbation modes that affect matter 
mixing in a spherical CCSN explosion. From the dynamical simulations of oxygen shell 


burning (Bazan & Arnett 1998 Meakin & Arnett 2007b), density perturbations of 5p/p ~ 
10% could be introduced at the boundaries of the convective shell burning. We assume 
that density perturbations are prominently introduced at the boundaries of shell burnings, 
and further assume perturbations with large amplitude (up to 50%) to obtain successful 
matter mixing. Then, we consider density perturbations at four critical radii (3 x 10 8 cm 
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at Si/C+O interface, 3 x 10 9 cm at oxygen shell burning position, 6 x 10 9 cm at C+O/He 
interface, and 5 x 10 10 cm at He/H interface), which roughly correspond to the boundaries 
of shell burnings in the progenitor. Density perturbation at each critical radius with the 
shape of 1 + e[2 rand(m#/7r) — 1] in the angular direction is given, where rand is the random 
number within the range between 0 and 1, and is a function of 9. We take m + 1 random 
numbers at sample points, 9 = 0 , 7r/m, 27T /m,..., (m — l)ir/m, tx. The random number for an 
arbitrary 9 is obtained by an interpolation with the values of adjacent sample points, and e 
is the amplitude of the perturbation. The amplitude has an exponential decay in the radial 
direction and can be presented as e = eo exp(—|r — r 0 |/r 0 ), where eo is the amplitude at each 
boundary layer and r 0 is the radius of the certain boundary. We call this case SCI. Thus, 
two parameters, eo and m, are involved in SCI. In order to understand the perturbation 
effect at each boundary on the matter mixing, we perform four more one-by-one boundary 
perturbation models, in which perturbations at the Si/C+O interface, the C+O/He interface, 
the He/H interface, and both C+O/He and He/H interfaces are considered, respectively. The 
RT instabilities at the H/He and He/C+O boundaries can be generated by the passage of 
the shock and/or reverse shock, and the condition of RT instability (Vp • Vp < 0) was given 
by Chevalier (1976). The growth rate of RT instabilities at H/He and He/C+O interfaces 
has been investigated (e.g., Ebisuzaki et al. 1989; Paper I) and such interfaces are unstable 
against RT instabilities after the passage of the shock waves. The four one-by-one boundary 
perturbation models are used to see which boundaries are critical for the matter mixing. In 
particular, we see the radial coherent structure in the perturbation pattern (Figure 2) that 
was obtained from large density perturbations imposed on those interfaces of the progenitor, 
which is important to the matter mixing of SN 1987A. 


Besides the possibility of large density perturbations at interfaces of different compo¬ 
sition layers, we propose that perturbations may also be introduced at any position in the 
progenitor. Although this global perturbation is not motivated by a known physical process, 
if the progenitor of SN 1987A was an LBV, such a perturbation might have existed. We 
propose the cases of SC2 and SC3, which are two different perturbation modes, to search 
the parameter space and explore the potential impacts. In the case of SC2, the density 
perturbation may be related to the density distribution of the progenitor star. We derive the 
density scale height H(r ), which can be presented as H{r ) = | dr/d lnp|. We assume that the 
radial length scale of the density perturbation / r (r) is proportional to the density scale height 
H{r). For simplicity, we introduce a scaling factor / of l r {r) = fH(r). If we adopt / = 5 
in a simulation model, we can obtain the radial coherent structure shown in Figure 2. A 
simulation model of / = 1 is also applied to compare with the corresponding model of / = 5. 
For a given radius, the angular length scale of the density perturbation is lg(r ) = rn/m, 
and we fix the number of m. We take random numbers at 2D sample points at intervals 
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of l r (r) and le{r). The perturbation at arbitrary (r, 9) is given by 1 + e 0 [2 rand(r, #) — 1], 
where the random number is obtained by an interpolation with the values of neighboring 
sample points. In the case of SC2, the perturbations are modulated not only in the angular 
direction but also in the radial direction. Thus, three parameters, eo, m and /, are involved 
in the case of SC2. 

Third, we assume that the length scales of the perturbation in both the radial and 
angular direction, l r (r) and lg(r), respectively, are proportional to the density scale liight 
H(r). In other words, we assume l r (r) = fH(r ) and lg(r) = fH(r). We take random 
numbers at 2D sample points at intervals of l r (r) and lg(r). The perturbation at arbitrary 
(r, 9) is given by 1 + eo[2 rand(r, 9) — 1], where the random number is obtained by an 
interpolation with the values of the neighboring sample points. We call this case SC3. For 
SC3, two parameters, eo, and /, are involvecQ 

The density scale hight, H(r), of the progenitor and the ratio between r and H(r) 
as a function of radius is illustrated in the upper panel of Figure 1. We can see that 
H(r) is roughly proportional to r and the ratio r/H(r) ranges between 2 and 8. The 
exception is the case for the regions around the stellar surface. In addition, H(r ) tends 
to be small around the composition interfaces of C+O/He and He/H. Moreover, we show 
the radial perturbations involved with the cases of SC2 and SC3 in the lower panel of 
Figure 1. We take the perturbation mode of SC3 (e 0 = 50% and / = 1) as an example, 
and see that the amplitude of the density perturbation in the radial direction is randomly 
distributed. Distributions of the 2D density perturbation clumps with amplitude eo = 50% 
for SCI, SC2, and SC3 are shown in Figure 2. For the case of SCI, perturbations are 
introduced around the critical radii, and angular perturbation is adjusted by the number of 
m. Perturbations extend from the chosen radii in both upward and downward directions. 
For SC2, perturbations are adjusted in angular direction and radial direction by m and 
/, respectively. For SC3, perturbations in both angular direction and radial direction are 
adjusted by the same number of /. Thus, from these perturbation patterns, we see that the 
coherent structure in the radial direction shown in some panels of Figure 2 plays a key role to 
reproduce high velocity 56 Ni clumps. In particular, we notice that the density perturbations 


1 Before we perform simulations of SCI, SC2, and SC3, one simple case, “Basicmodel”, is considered. In 

this case, we introduce density perturbations at the four critical radii in the progenitor star as mentioned 

in SCI. However, we set the amplitude of the density perturbation at each critical radius by the way of 

e = eoexp(—|r — ro\/h), where eo = 0-5, and h corresponds to the density scale height at each critical radius. 

From the function H{r ), we have h\ = 2.2 x 10 1 * * * * * * 8 cm, /i 2 = 1.2 x 10 9 cm, h 3 = 2.0 x 10 9 cm, and /14 = 2.0 x 10 10 

cm at each critical radius, respectively. Density perturbation at each critical radius in angular direction is 

as same as that given in SCI. In this simple case, we cannot reproduce high velocity 56 Ni clumps. Then, we 

turn to complicated cases of SCI, SC2, and SC3 in the spherical explosion to reproduce the matter mixing. 
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driven by the turbulent convection in stellar interiors show the coherent feature on large 
scales, and wave numbers have an order of the circumference at the outer boundary divided 


by the depth of the convection zone (Meakin & Arnett 2007a Arnett et al. 2009 Arnett 


& Meakin 2011a|b Viallet et ah 2013). The coherence shown in our perturbation patterns 


with m — 8 can capture a similar feature. 


2.3.2. A spherical Explosion 


The density perturbation modes mentioned in Section 2.3.1 are under the condition of 
spherical explosion. Because the observations of SN 1987A indicate an aspherical explosion, 
we can mimic some aspherical explosion cases that were adopted in Paper I. All of the 
following aspherical explosion models are under the case of SCI density perturbation mode. 

We take the initial radial velocity v r oc r[ 1 + acos(2d)]/(l + a) that was written as 
Equation (7) of Paper I. Through this function form, we can consider a bipolar explosion. In 
order to represent the bipolar explosion proposed by Nagataki et al. (1997, 1998), Nagataki 
(2000), and Paper I, we obtain the ratio of the radial velocity along the polar direction to 
that along the equatorial plane as u po i/v e q = (1 + a)/(l — a), where a is the parameter to 
determine the degree of asymmetry (a = 0 corresponds to a spherical explosion). In this 
paper, a is set to be 3/5 , i.e. v po \/v eq = 4. 

We also propose the asymmetric feature of the CCSN explosion to the equatorial plane 
as mentioned in Paper I (e.g., due to neutrino-driven explosion aided by SASI; Suwa et al. 
2010). An equatorially asymmetric explosion is mimicked by changing the normalization of 
the radial velocity across 6 = tt/2. We define the initial radial velocity at 6 = 0 and that at 
9 = 7 T as u up and t’d OW n , respectively. We adjust the normalization so that u up /t’down = 1.8. If 
we take u up /ud OW n = 2, as given by Paper I, the estimated neutron star kick velocity becomes 
larger than that typical of young pulsars. Therefore, we adopt a smaller value for u„ p /t'down 
in this paper. 


We can further combine the cases of bipolar explosion and equatorially asymmetric 
explosion together as a case of global asymmetric explosion. The combined global asymmetric 
explosion model may be used to examine SN 1987A explosion, as indicated by the asymmetric 
feature of [Fe II] line ( Haas et al.||1990 Spyromilio et al. 1 1990) and the elongated morphology 
(Kjaer et al. 2010) in SN 1987A. 
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3. Results 

Our 2D hydrodynamic simulation results on the matter mixing of CCSN explosions are 
presented in this section. We list those parameters adopted in our simulation models in 
Table 1. Some major simulation results are listed in Table 2. 


3.1. Spherical Explosion 

We first present our simulation results for SCI. The results of the case with the pertur¬ 
bation amplitude eo = 50% and rn — 8 are shown in Figure 3, whereas the results of the case 
with the perturbation amplitude eo = 50% and m — 20 are shown in Figure 5. We show den¬ 
sity distribution, distribution of the 56 Ni mass fraction, radial velocity of several elements, 
and line of sight velocity profiles of 56 Ni, respectively]^] We clearly see that prominent RT 
fingers are formed and distributions of 56 Ni have anisotropic structures. Although most 56 Ni 
is within the inner part of the contour, we see that the RT fingers of 56 Ni are beyond the 
radial velocity of 3,000 km s -1 . The matter mixing happens, as we can see in radial velocity 
profiles. The line of sight velocity of 56 Ni is extended beyond 3,000 km s _1 , and high velocity 
56 Ni clumps are reproduced. The perturbation modes of SCI show a prominent coherent 
structure in radial direction (Figure 2), which we believe is important to reproduce high 
velocity 56 Ni clumps. A different set of random numbers may induce different results. To 
test this, we use the model of SClp50m8 with a different set of random numbers and run the 
simulation again. We obtain the following results: total 56 Ni mass is 0.15 M 0 ; high velocity 
(> 3000 km s -1 ) 56 Ni is 21.7%; total 44 Ti mass is 2.1 x 10 -4 M 0 . Thus, our simulation 
results using this different set of random numbers are consistent with those of the original 
one. We plot the final density distribution, final distribution of 56 Ni mass fraction, radial 
velocity profile, and line of sight velocity profile of 56 Ni with different observational view 
angles in Figure 4. 

The perturbation mode in the angular direction is dominated by the parameter m. 
Figures 3 and 5 show that different values of parameter m provide different RT fingers: More 
RT fingers in the angular direction are developed if we take a larger value of m. However, 
as high velocity 56 Ni clumps are mainly determined by the strong coherent structure in the 
radial direction, a larger percentage number of high velocity 56 Ni can be obtained by the 
simulation models with a smaller number of m. 


2 In this paper, for simplicity, the mass distributions of elements as a function of radial velocity at the 
ends of simulation time are called radial velocity profiles, and mass distributions of elements as a function 
of line of sight velocity at the ends of simulation time are called as line of sight velocity profiles. 
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In order to investigate the effect of the perturbation amplitude on the matter mixing 
of SN 1987A, simulations with 10% of density perturbation are also performed in the case 
of SCI. We show the results using parameter m = 20 as an example in Figure 6. Density 
distribution, distribution of 56 Ni mass fraction, radial velocity profiles of several elements, 
and line of sight velocity profiles of 56 Ni are given, respectively. Although RT instability is 
developed at the last stage, 56 Ni has roughly spherical distribution. We see that line of sight 
velocity of 56 Ni is limited within 2,000 km s -1 , and the mixing of the innermost elements 
(e.g., 56 Ni and 44 Ti) into high velocity (> 3,000 km s _1 ) regions has failed. We conclude 
that a density perturbation as small as 10% cannot reproduce the high velocity 56 Ni clumps 
of SN 1987A. 

It is interesting to examine which boundary layer where density perturbation occurs is 
important for the matter mixing. We perform three simulation models in which perturbations 
are introduced at a single composition (Si/C+O, C+O/He, He/H) interface. One simulation 
model with perturbations at both C+O/He and He/H interfaces is also performed. The 
parameters of m = 20 and eo = 50% are set for each simulation model. We do not obtain 
high velocity 56 Ni when the density perturbation is introduced at the boundary layer between 
silicon and oxygen. If density perturbation is introduced at the C+O/He interface or He/H 
interface, high velocity 56 Ni clumps appear. However, the mass of high velocity 56 Ni is 
only about 0.7% of the total 56 Ni produced by the CCSN explosion, which is less than the 
observational value (> 4%, Haas et al. 1990). If we take perturbations at both C+O/He and 
He/H interfaces, high velocity clumps of 56 Ni are successfully obtained. These results indicate 
that most effective density perturbations for the matter mixing occur at both C+O/He and 
He/H interfaces, which is consistent with the examination shown in Figure 4 of Paper I—the 
growth factors of seed perturbation around the interfaces of C+O/He and He/H in the 
progenitor model are prominent. Figure 7 shows results obtained from the model in which 
the perturbations at C+O/He and He/H interfaces are introduced. 

We then present our simulation results for SC2. The results of the case with perturbation 
amplitude eo = 50% and m — 8 are shown in Figure 8, while the results of the case with 
perturbation amplitude eo = 50% and m = 20 are shown in Figure 9. For both cases, we 
take parameter / = 5. We see that RT fingers are developed and the matter mixing works: 
We obtain high velocity 56 Ni clumps. Cases with a small perturbation amplitude eo = 10% 
are also presented, and the results obtained from rn — 8 and m = 20 are given. We show the 
results from the case of e 0 = 50%, m — 20, and / = 5 as an example in Figure 10. Although 
we see some developed RT fingers, we confirm that the SC2 model with such small density 
perturbation cannot reproduce the matter mixing. 

We further investigate the effect of radial coherent structure on the matter mixing 
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for SC2. This effect is dominated by the parameter /. A small number of / indicates a 
small length scale of density perturbation. Figure 11 shows the results of the simulation of 
eo = 50%, m = 20, and / = 1 as an example. We cannot reproduce the matter mixing 
by this simulation model. As we can see from the density perturbation patterns in Figure 
2, the perturbation mode of / = 5 provides the prominent coherent structure in the radial 
direction, and the matter mixing is obtained. The perturbation mode of / = 1 provides a 
far less coherent structure in the radial direction, so the matter mixing cannot be obtained. 

We present our simulation results for SC3 at last. The density perturbation with 50% 
amplitude is considered in two simulation models, and / = 1 and / = 5 are taken, respec¬ 
tively, in each one. However, we do not find high velocity 56 Ni clumps in these simulation 
results. Figure 12 shows the results from the simulation model with eo = 50% and f — 1. As 
we see from the density perturbation patterns in Figure 2, the perturbation modes of SC3 
do not provide any clear coherent structure in the radial direction. Thus, the matter mixing 
cannot be obtained. 

We summarize our main results of the spherical explosion models: (1) The density 
perturbation amplitude is the most important parameter and RT instability dominates the 
matter mixing. RT fingers can be developed and the matter mixing can be reproduced by 
some simulation models with 50% of density perturbation. While 10% of density pertur¬ 
bation introduced in simulation models can generate RT fingers, but it is not enough to 
reproduce the matter mixing. (2) The synthesized 56 Ni and 44 Ti in the stellar interior can 
be effectively dug out because the radial coherent structures shown in perturbation patterns 
make the successful matter mixing in the CCSN explosion. (3) The matter mixing cannot be 
reproduced by the simulation model with the density perturbation introduced at the Si/C+O 
interface, but density perturbations at both the C+O/He and He/H interfaces can effectively 
achieve high velocity 56 Ni clumps. In other words, prominent density perturbations at large 
scales are successful for the matter mixing of a CCSN explosion. (4) The distribution of 
56 Ni is roughly spherical in the simulation models where the matter mixing cannot be re¬ 
produced. If the simulation models are successful for the matter mixing, the distribution 
of 56 Ni is expected to be strongly anisotropic with prominent RT fingers. The properties of 
the heavy elements distribution shown in our simulation results indicate that SN 1987A may 
not favor spherical explosion models. (5) The distribution of 44 Ti follows that of 56 Ni. (6) 
Although we initially assume a spherical explosion in our simulation models, large density 
perturbations introduced in the progenitor may break down this spherical property. The set 
of initial conditions is important because it affects the evolution of the CCSN explosion. 

We focus on 56 Ni velocity profiles in these spherical explosion models, and obtain high 
velocity 56 Ni clumps. We think that the matter mixing of SN 1987A can be qualitatively 
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reproduced by these models. We also further discuss two issues. First, the total mass of 
synthesized 56 Ni in our spherical simulation results is about two to three times of that of 
observational determinations (Shigeyama et al. 1988). Second, the explosion energy pro¬ 
duced by our spherical simulations is about two times that of the observational fittings by 
Shigeyama et al. (1988[) and Shigeyama & Nomoto (1990), although the explosion energy 


has a range of (1 — 2) x 10 51 erg (Handy et al. 


2014). We note that Kifonidis et al. (2006) 


produced the matter mixing in the CCSN with an explosion energy of 2 x 10 51 erg, and 
Hammer et al. (2010) provided an explosion energy of 1.8 x 10 51 erg in their 2D models. 


In order to examine the matter mixing of CCSN spherical explosions affected by the input 
energy, we reduce the input energy from 2.5 x 10 51 to 1.7 x 10 51 erg and re-perform one 
simulation model for SCI. Those results are shown in Figure 13. The explosion energy of 
1.08 x 10 51 erg calculated from this simulation model is consistent with the observational 
one given by Shigeyama et al. (1988). But the total synthesized 56 Ni mass is still larger than 
the observational one, and we do not obtain enough high velocity 56 Ni clumps. Therefore, 
it is necessary to perform simulation models of aspherical CCSN explosion to study matter 
mixing. 


3.2. Aspherical Explosion 

We present the results of aspherical explosion models in this subsection. First, we show 
the results of a bipolar model. Density distribution, distribution of the 56 Ni mass fraction, 
radial velocity of several elements, and line of sight velocity profiles of 56 Ni are shown in 
Figure 14. Due to this jet-like CCSN explosion, newly synthesized heavy elements in the 
stellar interior are pushed outside, mainly along the polar direction. Distribution of 56 Ni is 
roughly along the polar regions, and we see from the radial velocity profiles that the matter 
mixing works. We also successfully obtain high velocity 56 Ni clumps. 

The results of the equatorially asymmetric explosion model are shown in Figure 15. 
Due to the asymmetric explosion across the equatorial plane, the density distribution and 
the distribution of 56 Ni are both equatorially asymmetric. The matter mixing is effective 
and high velocity 56 Ni clumps are obtained. 

We also obtain the results of the combination model from the bipolar explosion and 
equatorial asymmetric explosion. In this global asymmetric case, the matter mixing is suc¬ 
cessfully produced. However, we overestimate masses of high velocity 56 Ni and 44 Ti when 
we introduce 50% of the density perturbation in the progenitor. In order to reproduce the 
proper high velocity 56 Ni clumps, we reduce the density perturbation amplitude from 50% 
to 25%. The resulting mass and velocity profiles of 56 Ni from the simulation model with the 
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reduced perturbation amplitude are consistent with the observational ones of SN 1987A. We 
also perform one simulation model of a global asymmetric case with 10% density perturba¬ 
tion amplitude. We confirm that the obtained high velocity 56 Ni clumps are not enough if 
only 10% of the density perturbation is introduced in the progenitor. 


We summarize our main results for the matter mixing in an aspherical CCSN explosion: 

(1) Either bipolar explosion or equatorially asymmetric explosion can effectively reproduce 
high velocity 56 Ni clumps if we introduce 50% of the density perturbation in the progenitor. 

(2) We cannot obtain enough high velocity 56 Ni clumps from the simulation model with 
10% of the density perturbation; at least 25% of the density perturbation is necessary to 
reproduce the observed matter mixing and high velocity 56 Ni clumps. (3) The distribution 
of 44 Ti follows that of 56 Ni. (4) The total 56 Ni mass and the mass of high velocity 56 Ni 
produced in aspherical explosion models roughly fit the observational values. However, the 


explosion energy is still overestimated compared to the number given by Shigeyama et al. 


(1988) and Shigeyama & Nomoto (1990). 


In order to examine the matter mixing of the CCSN aspherical explosion affected by 
the input energy, we can reduce the input energy and re-perform simulation models. For 


example, we obtain our simulation results if the input energy is 2.0 x 10 51 erg: the total 
56 Ni mass is O.O7M 0 , high velocity (> 3,000 km s _1 ) 56 Ni mass is 6.5% of the total 56 Ni 
mass, and explosion energy is 1.43 x 10 51 erg. The density distribution, distribution of 56 Ni 
mass fraction, radial velocity profiles of several elements, and line of sight velocity profiles of 
56 Ni are shown in Figure 16. These results are consistent with the observational ones of SN 
1987A. We regard this model (Asph+SClp25m20** listed in Tables 1 and 2) as one of the 
most favorable simulation models in our work to explain the observational phenomena of the 
matter mixing in SN 1987A. The bipolar feature inferred from rotational CCSN explosion 
has been also explored by theoretical models (Nakamura et al. 2014). 


4. Discussion 


We introduce large density perturbations in the CCSN progenitor to investigate the 
matter mixing of the CCSN explosion. The consequent question is: What is the mechanism 
for the large density perturbation that occurs in stellar interiors? MLT is usually considered 
in the convection zones of the stellar interior (Kippenhahn & Weigert 1990). However, MLT 
is a local and spherically symmetric treatment in stellar models. Non-local and anisotropic 


convective mixing theories have been developed ( 

Xiong 

1977 

Grossman et al. 

1993; 

Xiong 

et ah 1997 

Canuto & Dubovikov 1998; 

Deng et al. 2006; 

Li & Yang 2007 

Hotta et al. 2014), 


but those treatments focus on the evolution of main-sequences stars. In addition, it seems 
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that the related overshooting in solar models and post-main-sequence star models (Tian et 


al. 2009 Li 2012; Zhang & Li 2012a|b ) cannot strongly affect the matter mixing of CCSN 


explosion because the overshooting model is applied for the globally stable state of stars, 
ffowever, the large density perturbation proposed in this paper is one requirement for the 


successful matter mixing of the CCSN explosion. Meakin & Arnett (2006) noticed that the 


spherically symmetric features in stellar interior can be broken as density perturbation is 
induced by turbulence within stellar convection zones. This kind of perturbation is related 
to carbon and oxygen shell burning. Analytic explorations and 3D simulations of turbulent 
convection have been performed, where turbulent entrainment and turbulent kinetic energy 


have been considered (Meakin & Arnett 2007a Arnett et al. 2009). About 10% of density 


perturbation is obtained in convection zones as carbon and oxygen shell burning (Meakin 


& Arnett 2007b). Recently, Arnett & Meakin ( 2011a|b ) found that violent dynamic erup¬ 
tions occur in convective zones with simultaneous C, O, and Si shell burning of a CCSN 
progenitor star via 3D simulations. A realistic stellar model should provide aspherical prop¬ 


erties from hydrodynamic instabilities and turbulence for CCSN explosion (Smith & Arnett 


2014). Moreover, it has been found recently that the evolution of matter mixing is different 


between an BSG progenitor model and an RSG progenitor model (Wongwathanarat et al. 


2015). In our 2D hydrodynamic simulation models, although 10% of density perturbation 


in the progenitor can produce RT fingers, 50% the of density perturbation is required in 
order to reproduce high velocity 56 Ni clumps. Our results provide additional clues to further 
analyze the turbulent convection of the CCSN explosion. 

We have performed 2D hydrodynamic simulations to research of CCSN matter mixing 


in this work. Kane et al. (2000) pointed out that RT instability grows faster in 3D simu¬ 


lation models, andlHammer et al. (2010) reached a similar conclusion after comparing 2D 


and 3Dsimulation results. As Hammer et al. (2010) explained, heavy element clumps have 


toroidal topology in 2D simulations and bubble-sphere topology in 3D simulations, so these 
heavy element clumps experience less drag force in 3D simulations than in 2D ones. Thus, 
heavy element clumps can retain higher velocities in 3D models. If these physical processes 
mentioned above can really happen in CCSN explosions, we expect that less density per¬ 
turbation should be introduced in 3D simulation models to reproduce high velocity 56 Ni 


clumps. However, it seems that Joggerst et al. (2010b) did not find qualitative differences 


of matter mixing efficiency between 2D and 3D simulations. One possible reason is that the 
interaction of the small RT fingers may reduce the final growth of RT instability, even if RT 
instability grows faster at the beginning of a CCSN explosion in 3D simulations. We can 
use 3D simulations with different density perturbation amplitudes to further investigate this 
CCSN matter mixing issue in the future. 


We notice that RT fingers do not show a prominent feature at the early stage of CCSN 
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explosion. Couch & Ott (2013) proposed velocity perturbations in an angular direction at the 
Si/O interface of a progenitor. Their 3D simulations showed that this kind of perturbation 
can affect CCSN explosion. If we consider the neutrino-driven convection, the CCSN shock 
will be asymmetric around Fe/Si or Si/O boundaries. In fact, Kifonidis et al. (2006) claimed 
that the initial global deformation of the shock can trigger strong instability at the He/H 
boundary, which makes the global anisotropy of the inner ejecta. We find that large density 
perturbations placed at the Si/O interface may only provide and/or enhance aspherical 
features at the beginning of the CCSN explosion, but this does not give any direct effects 
on the matter mixing at the later time of the CCSN explosion. This discrepancy may be 
due to the application of different progenitors in different hydrodynamic simulation models. 
In our simulations, RT fingers are clearly identified after the shock propagates though the 
C+O/He interface. Giant RT fingers are developed quickly after the forward shock passes 
the He/H interface, and finally dominate at the late stage of CCSN explosion. During this 
later explosion process, we see a few small RT fingers merging into a large one, which was 
also mentioned by other previous works (Arnett et al.||1989 Kifonidis et al. [2003 Hammer et 


al. 2010 Ellinger et al. 2012). Moreover, when we successfully reproduce high velocity 56 Ni 


clumps, we clearly see the huge and prominent RT fingers with the length scale of about 
2 x 10 12 cm in these simulation models. However, the RT fingers with a length scale smaller 
than 1 x 10 12 cm shown in some simulation models are short and weak, and we cannot 
reproduce high velocity 56 Ni clumps by these models. 


The obtained total 56 Ni mass from spherical explosion models in this paper tends to be 
overestimated compared with that of SN 1987A. There are two possible reasons: First, if we 
reduce the input energy in spherical explosion models, the total 56 Ni mass can be reduced; 
second, a limited nuclear reaction network is included in the FLASH code and some isotopes 
of 56 Ni are neglected. Although a detailed discussion on the production of isotopes goes far 
beyond the scope of this paper, we attempt to examine the effect of neglecting isotopes on the 
production of 56 Ni by using a large nuclear reaction network. We perform a ID hydrodynamic 
simulation of CCSN explosion without any perturbation using the same progenitor model 
mentioned in Section 2. The input energy is set to be 2.5 x 10 51 erg. Then, explosive 
nucleosynthesis is calculated as a post-process using the large nuclear reaction network, 
including 464 nuclei. As a result, the following values are obtained: M( 56 Ni) = 0.14 M 0 , 
M( 5 'Ni) = 0.01 M 0 , and M( 58 Ni) = 0.05 M & . On the other hand, M( 56 Ni) obtained by ID 
FLASH code is 0.21 M 0 . This indicates that the nuclear reaction network in the FLASH 
code may overestimate the mass of 56 Ni in a CCSN explosion. We present the detailed 
calculation process using the large nuclear reaction network in the Appendix. 


We also notice that the masses of high velocity 56 Ni and 44 Ti clumps could be affected 
by the dimensionality in this paper. It is declared that 2D axisymmetric hydrodynamic 
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simulations have numerical problems due to the computational grid-geometry. Kifonidis et 
al. ( 2006[ ) and Gawryszczak et al. (2010) mentioned that RT fingers may grow faster along 
the polar direction due to the discretization errors and the non-penetrating feature around 
the polar axis. Thus, the artificial axial jets are shown in axisymmetric simulation results. 
We have also observed these kinds of features in our simulation contours. We speculate 
that some extra high velocity clumps of 56 Ni and 44 Ti along the polar direction may be 
artificially produced by these axial jets. The productions of total 56 Ni and 4i Ti may be 
also affected by this geometric effect. In order to examine the axial jet effect, we perform a 
simple test. We consider the computational domain with angular direction from 15 degree 
to 165 degree and calculate total 56 Ni mass and high velocity 56 Ni clumps. We call this 
case “excludel5”. For comparison, we also consider the computational domain with angular 
direction from 30 degree to 150 degree and calculate the same quantities again. We call this 
case “exclude30”. Two simulation models-the spherical explosion model “SClp50m20” and 
the bipolar explosion model “Bipo+Clp50m20”-are chosen. Line of sight velocity profiles and 
radial velocity profiles are shown in Figure 17. We see that the axial jet effect may produce 
artificial high velocity 56 Ni clumps. Bipolar aspherical explosion models may suffer from this 
effect more seriously than spherical explosion models. The related quantities are as follows: 
for the case of “excludel5” and “SClp50m20”, total 56 Ni mass is 0.16 M 0 and high velocity 
(> 3, 000 km s 4 ) 56 Ni clumps have 5.5% of the total 56 Ni mass; for the case of “exclude30” 
and “SClp50m20”, total 56 Ni mass is 0.15 M 0 and high velocity (> 3,000 km s~ 4 ) 56 Ni 
clumps have 5.1% of the total 56 Ni mass; for the case of “excludel5” and “Bipo+Clp50m20”, 
total 56 Ni mass is 0.08 M 0 and high velocity (> 3, 000 km s^ 1 ) 56 Ni clumps have 8.5% of the 
total 56 Ni mass; for the case of “exclude30” and “Bipo+Clp50m20”, total 56 Ni mass is 0.06 
Mq and high velocity (> 3,000 km s _1 ) 56 Ni clumps have 0.9% of the total 56 Ni mass. As a 
summary, the axial jet effect may only reduce a few percent of high velocity 56 Ni clumps in 
a spherical explosion model, but it may reduce at least 20% of high velocity 56 Ni clumps in 
a bipolar explosion model. 3D simulations are necessary for further quantitative analysis. 

The observed line profiles of [Fe II] in SN 1987A are not symmetric and their peak is 


redshifted (Haas et al. 1990). The inner ring of SN 1987A is inclined at about 45° to the 


sky: the northern part of the ring is close to observers and the southern part of the ring is 


away from observers (Tziamtzis et al. 2011). The inner ejecta of SN 1987A is elongated and 


it is roughly in the same plane of the inner ring (Kjaer et al. 2010). We calculate the line 
of sight velocity profiles of 56 Ni to different observational view angles. We define positive 
velocity and negative velocity in the line of sight velocity profiles as the redshifted velocity 
and blueshifted velocity, respectively. From our analysis of density contours and line of sight 
velocity profiles of 56 Ni, it looks like that the profile with observational view angle of 45° 
or 135° is morphologically close to the observational one. In particular, the peak shift is 
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dominated by the ratio of n up and Wdown when we consider equatorially asymmetric explosion 
models. A large value of the ratio leads to a strong peak shift. On the other hand, the 
density perturbation also affects the peak shift in equatorially asymmetric explosion models. 
We find that the large density perturbation tends to weaken the peak shift. Some fine 
tuning in simulation models is required to reproduce exactly the observational redshifted 
peak of line profile, although this is not the main scope in this paper. The results of our 2D 
simulations are helpful to understand the morphology of SN 1987A, but it is obvious that 
3D simulations, which we plan to perform in the future, are more suitable to explain the 
reconstructed 3D features of SN 1987A flKjaer et~ah 2010; Larsson et al. 2013). 

Direct detections of emission lines from the decay of 44 Ti in two young supernova rem¬ 
nants have been reported. The mass of 44 Ti in SN 1987A was derived as (3.1±0.8) x 10 ~ 4 M 0 


Grebenev et al. 

2012 

) and the mass of 44 Ti in Cassiopeia A is about 1 . 6+03 x 10 4 M 0 

'Grefenstette et al. 

2014 

). We also obtain the mass of 44 Ti in our simulation models. Our 


values obtained from the spherical explosion models are roughly comparable to the observa¬ 
tional one, while the values obtained from the aspherical explosion models are larger than 
the observational one. Moreover, the ratios between the 56 Ni mass and 44 Ti mass in spherical 
explosion models are larger than those of aspherical models. If we overestimate the mass of 
56 Ni in the spherical explosion models, the process of an a-rich freeze-out may enhance the 
production of 44 Ti ( Woosley fc Hoffman|l991 Nagataki|2000 ). The so-called mass-cut, which 
arbitrarily determines the location of the ejecta from the compact remnant may affect the 
mass of 4 i Ti. However, it seems in Paper I that the differences of velocity profiles between the 
results with mass-cut and those without mass-cut are not significant. We also note that the 
overestimation of both 56 Ni and 44 Ti may be induced by the small nuclear reaction network 
in the FLASH code. Compared to the productions from spherical explosion models, it seems 
that those obtained from aspherical explosion models are more consistent with observational 
results. We confirm the speculation of Paper I that 44 Ti may be an indicator of explosion 
asphericity. However, we also realize that the production of 4i Ti and 56 Ni in the CCSN 
explosions is related to not only temperature and density, but also the electron fraction and 
detailed thermodynamic trajectories (Magkotsios et al. 2010). Postprocess calculations of 
the nucleosynthesis with a large nuclear reaction network are required. These complicated 
physical processes make some additional uncertainties to the ratio between 44 Ti and 56 Ni. 


We can calculate radial velocities not only for the elements of H, 56 Ni, and 4 i Ti, but 
also for the elements of 4 He, 12 C, 16 0, and 28 Si. The radial velocities of the elements in 
each simulation model are summarized in Table 3. If the simulation models of the spherical 
explosion with large density perturbation (50% amplitude) can produce high velocity 56 Ni 
and 44 Ti, the radial velocity profiles of 12 C, 16 0, and 28 Si also extend to a large velocity 
part. On the other hand, the simulation models of a spherical explosion with small density 
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perturbation (10% amplitude) cannot produce high velocity 56 Ni and 14 Ti, and it is also 
difficult to obtain more high velocity heavy elements of 12 C, 16 0, and 28 Si. In order to 
examine the effect of large density perturbation on the distributions of 4 He, 12 C, 16 0, and 
28 Si, we plot the line of sight velocity profiles of these elements in Figure 18 (view angle 
90 degree) and Figure 19 (view angle 45 degree). Three simulation models (SClp50m20, 
SClpl0m20, and Asp+Clp50m20(*)) are considered. We see that 12 C, 16 0, and 28 Si are 
extended into the high velocity region in the spherical explosion model with 50% of density 
perturbation. Aspherical explosion with the view angle of 90 degree takes effects on the 
distribution of 12 C, 16 0, and 28 Si, so that these elements are concentrated in the center of 


line of sight velocity profiles. From the observational point of view, Larsson et al. (2013) 


have obtained [Si I] + [Fe II] line profiles of SN 1987A. The [Si I] + [Fe II] emission is close to 
the plane of the ring and the profile is north-south asymmetry that is partially due to dust 
absorption. The emission is concentrated in the velocity interval of 1500 — 3000 km s -1 . 
Tanaka et al. (2012) have measured the polarization features of [He I], [O I], [Ca II], and [Fe 


II] in type Ib/Ic SN 2009jf and SN 2009mi. In principle, this kind of spectral observation 
provides the opportunity for us to extract the velocity profiles of these elements. We expect 
that those velocity profiles calculated from our simulation models can be further examined by 
future observations and consequently provide additional hints to explore clumpy distributions 
and multi-dimensional geometries of CCSN explosions. 


5. Summary 

We perform 2D hydrodynamic simulations of a CCSN to explain the high velocity 56 Ni 
clumps of SN 1987A. The density perturbation with 50% amplitude is introduced in the 
progenitor model, such that the matter mixing can be reproduced in spherical explosion 
models. We take this large density perturbation to be one possible way to explain the final 
aspherical feature of a CCSN explosion indicated by the distributions of 56 Ni and 4i Ti. On 
the other hand, simulation models with 10% of density perturbation cannot reproduce high 
velocity 56 Ni clumps. The RT instability is efficient at both the C+O/He interface and He/H 
interface to reproduce high velocity 56 Ni clumps. The radial coherent structures shown in 
perturbation patterns are important to produce the matter mixing. 

When we consider the case of aspherical explosion, at least 25% of density perturbation 
introduced in the progenitor is necessary to reproduce the high velocity 56 Ni clumps of SN 
1987A. Either bipolar explosion or equatorially asymmetric explosion is effective for the 
explanation of observational matter mixing properties. 

Explosion energy, total 56 Ni mass, total 44 Ti mass, and high velocity (> 3,000 km s _1 ) 
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56 Ni mass have been obtained in each CCSN simulation model. The distribution of 44 Ti 
follows that of 56 Ni in all simulations models. After the comparison of our simulation results 
to the observations of SN 1987A, we note that one of the most favorable simulation models 
in this work is an aspherical one with the combination of bipolar explosion and equatori- 
ally asymmetric explosion, in which 25% of density perturbation is introduced at different 
composition interfaces of the CCSN progenitor. 


As noted in Section 4, we have performed 2D simulations for the research of the matter 
mixing. However, 3D simulations may result in different matter mixing compared with 2D 
ones (Hammer et ah 2010). If we perform 3D simulations using the same progenitor, we 
expect that less density perturbations will be introduced to make the same matter mixing 
of SN 1987A. To obtain a more robust conclusion, we will perform 3D simulations for the 
matter mixing of SN 1987A in the future. 


We are grateful to Dr. A., Wongwathanarat and T. Takiwaki for their useful discussion. 
The software used in this work was in part developed by the DOE NNSA-ASC OASCR 
Flash Center at the University of Chicago. Numerical computation in this work was carried 
out on SR16000 at the Yukawa Institute Computer Facility. Mao, J. is supported by the 
Hundred-Talent Program (Chinese Academy of Sciences). This work is also supported by the 
Ministry of Education, Culture, Sports, Science and Technology (MEXT: No. 23105709 and 
No. 26105521) and the Japan Society for the Promotion of Science (JSPS: No. 19104006, 
No. 23340069, No. 2503786, No. 25610056, No. 26287056, No. 26800159, No. 26800141 
and No. 24.02022). 


A. Production of 56 Ni and 44 Ti with a Large Nuclear Reaction Network 


CCSN explosive nucleosynthesis has been comprehensively studied (Woosley & Arnett 


1973 Hashimoto 1995) and applied for SN 1987A (Hashimoto et al. 1989 Woosley & Hoffman 


1991). In this paper, we perform 2D hydrodynamic simulations using the FLASH code to 


study the matter mixing of CCSN explosion. As noted in Nagataki et ah (1997), some 
isotopes of 56 Ni should be included in the explosive nucleosynthesis. In order to clarify the 
effect of the small nuclear reaction network in the FLASH code to the production of 56 Ni 
and 44 Ti, we calculate the explosive nucleosynthesis using a large nuclear reaction network 
and a fiducial pure spherical explosion model. 


We first perform a ID hydrodynamic simulation of CCSN explosion using the same 
progenitor model and methods. The input energy of 2.5 x 10 51 erg is adopted for a fiducial 
case. To calculate explosive nucleosynthesis as a post-process, we adopt the so-called tracer 
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particle method (Nagataki et al. 1997): 2,000 particles are distributed on the computational 
domain of the hydrodynamic simulation (1.36 x 10 8 cm < r < 6 x 10 9 cm, which covers up to 
the oxygen-rich layer), and mass and composition are assigned to each particle to reproduce 
the mass distribution and the composition of the progenitor model. Then, by using the 
results from the time evolution of the velocity held in Eulerian grids, particles are passively 
moved to obtain the time evolution of density and temperature. Using the initial composition 
and the time evolution of density and temperature, explosive nucleosynthesis is calculated 
for each ejected particle. We adopt a large nuclear reaction network including 464 nuclei (up 
to 94 Kr), which was used in Ono et al. (2009). It should be noted that the innermost region 
(r = 1.36 x 10 8 cm) is slightly neutron-rich (electron fraction Y e ~ 0.485) and the production 
of neutron-rich isotopes, 57 Ni and 58 Ni, is expected. As a result, we obtain the following 
values: M( 56 Ni) = 1.4 x 10~ 4 M 0 , M( 57 Ni) = 1.0 x 10~ 2 M 0 , M( 58 Ni) = 5.4 x 10~ 2 M 0 , and 
M( 44 Ti) = 1.0 x 10~ 4 M 0 . The masses of 56 Ni and 44 Ti estimated by the FLASH code are 
2.1 x 10 _1 M q and 3.0 x 10~ 4 M 0 respectively. Thus, the masses of 56 Ni and 44 Ti calculated 
by the FLASH code are overestimated by a factor of 0.7 and a factor of 3, respectively, 
compared with those calculated by the large nuclear reaction network. Because the amount 
of 44 Ti is much less than that of 56 Ni, the overestimation of 44 Ti due to neglecting other 
elements is more severe than that of 56 Ni. 
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Fig. 1.— Upper panel: density scale height H(r) and r/H(r) of the progenitor model as a function of 
radius; lower panel: profiles of density and density perturbation in the progenitor model. Distributions of 
original density, perturbative density, and perturbation amplitude are denoted by red, blue, and cyan lines, 
respectively. Perturbation amplitude is randomly distributed in the radial direction and the maximum value 
is 50%. We take this perturbation mode from the case of SC3 (e 0 = 0.5 and / = 1). 
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Kig. 2.- Distributions of density perturbation in the r-9 plane. The ranges of r and 9 are 1.4 x 10 8 cm < r < 3.4 x 10 12 cm 

and 0 < 9 < 7r, respectively, and the color bar in the top of each panel denotes the amplitude range of perturbation, upper 
left: perturbation mode of SCI (m=8); Upper right: perturbation mode of SCI (m=20); middle left: perturbation mode of 
SC2 (m = 20 and / = 1); middle right: perturbation mode of SC2 (m = 20 and / = 5); lower left: perturbation mode of SC3 
(m = 20 and / = 1); lower right: perturbation mode of SC3 (m = 20 and / = 5). We note four critical radii: 3 x 10 8 cm at 
Si/C+O interface, 3 X 10 9 cm at oxygen shell burning position, 6 x 10 9 cm at C+O/He interface, and 5 X 10 10 cm at He/H 
interface. Only the models with coherent-structured perturbation in the radial direction can reproduce high velocity clumps of 
56 Ni (see upper left, upper right, and middle right panels). 
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Fig. 3.— Results obtained from the model SClp50m8. Upper left: final density distribution; 
upper right: final distribution of 56 Ni mass fraction. The white and black lines indicate the 
contours of radial velocities of 3,000 km s -1 and 4,000 km s^ 1 , respectively; Lower left: radial 
velocity profile; lower right: line of sight velocity profile of ,a6 Ni with different observational 
view angles. For each lower panel, the velocity bin is 100 km s -1 , numbers of A M/M 
are normalized by their maximum value, where AM is the mass in the velocity range of 
v ~ v + An; and M is the total mass. 
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Fig. 4.— Results obtained from the model SClp50m8 with a different set of random numbers. 
Panel descriptions are the same as those in Figure 3. 
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Fig. 5.— Results obtained from the model SClp50m20. Panel descriptions are the same as 
those in Figure 3. 
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Fig. 6.— Results obtained from the model SClpl0m20. Panel descriptions are the same as 
those in Figure 3. 
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Fig. 7.— Results obtained from the model SCIC+O/He/H. Panel descriptions are the same 
as those in Figure 3. 
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Fig. 8.— Results obtained from the model SC2p50m8f5. Panel descriptions are the same as 
those in Figure 3. 
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Fig. 9.— Results obtained from the model SC2p50m20f5. Panel descriptions are the same 
as those in Figure 3. 
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Fig. 10.— Results obtained from the model SC2pl0m20f5. Panel descriptions are the same 
as those in Figure 3. 
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Fig. 11.— Results obtained from the model SC2p50m20fl. Panel descriptions are the same 
as those in Figure 3. 
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Fig. 12.— Results obtained from the model SC3p50fl. Panel descriptions are the same as 
those in Figure 3. 
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Fig. 13.— Results obtained from the model SClp50m20(*). Panel descriptions are the same 
as those in Figure 3. 
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Fig. 14.— Results obtained from the model Bipo+SClp50m20. Panel descriptions are the 
same as those in Figure 3. 
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Fig. 15.— Results obtained from the model AsyEqua+SClp50m20. Panel descriptions are 
the same as those in Figure 3. 
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Fig. 16.— Results obtained from the model Asph+SClp50m20(**). Panel descriptions are 
the same as those in Figure 3. 
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Fig. 17.— Line of sight velocity profiles and radial velocity profiles of 56 Ni for the simple 
test of the axis-jet effect. The labeled value in the y-axis of each panel is A M/M (without 
normalization), where AM is the mass in the velocity range of v ~ v + An, and M is the 
total mass. The view angle is 45 degree and the velocity bin is 100 km s -1 . Upper left: 
Line of sight velocity profiles of the model “SClp50m20”; upper right: line of sight velocity 
profiles of the model “Bipo+Clp50m20”; lower left: radial velocity profiles of the model 
“SClp50m20”; lower right: radial velocity profiles of the model “Bipo+Clp50m20”. In each 
panel, “0—180”, “15—165”, and “30—150” indicate the computational domain within the 
angular ranges of 0—180 degree, 15—165 degree, and 30—150 degree, respectively. 
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Fig. 18.— Line of sight velocity profiles of four elements. The view angle is 90 degree and 
the velocity bin is 100 km s _1 . Upper left: line of sight velocity profiles of 4 He; upper right: 
line of sight velocity profiles of 12 C; lower left: line of sight velocity profiles of 16 O; lower 
right: line of sight velocity profiles of 28 Si. The labeled value in the y-axis of each panel is 
A M/M (without normalization), where AM is the mass in the velocity range of v ~ v + Av, 
and M is the total mass. Three models are considered in each panel: SClp50m20 (spherical 
explosion with 50% of density perturbation), SClpl0m20 (spherical explosion with 10% of 
density perturbation), and Asp+Clp50m20(*) (aspherical explosion with 50% of density 
perturbation). 
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Fig. 19.— Line of sight velocity profiles of four elements. The view angle is 45 degree and 
the velocity bin is 100 km s _1 . Upper left: line of sight velocity profiles of 4 He; upper right: 
line of sight velocity profiles of 12 C; lower left: line of sight velocity profiles of 16 O; lower 
right: line of sight velocity profiles of 28 Si. The labeled value in the y-axis of each panel is 
A M/M (without normalization), where AM is the mass in the velocity range of v ~ v + Av, 
and M is the total mass. Three models are considered in each panel: SClp50m20 (spherical 
explosion with 50% of density perturbation), SClpl0m20 (spherical explosion with 10% of 
density perturbation), and Asp+Clp50m20(*) (aspherical explosion with 50% of density 
perturbation). 

























Table 1: Parameters in simulation models 


Model/Parameter 

e 0 

m 

/ 

^pol / ^eq 

^up/^down 

Perturbation Position 

Fiducialmodel 

0.0 

0 

- 

- 

- 

no perturbation 

Basicmodel 

0.5 

20 

- 

- 

- 

4 critical radii 

SClp50m8 

0.5 

8 

- 

- 

- 

4 critical radii 

SClp50m8(O) 

0.5 

8 

- 

- 

- 

4 critical radii 

SClp50m20 

0.5 

20 

- 

- 

- 

4 critical radii 

SClp50m20(*) 

0.5 

20 

- 

- 

- 

4 critical radii 

SClpl0m8 

0.1 

8 

- 

- 

- 

4 critical radii 

SClpl0m20 

0.1 

20 

- 

- 

- 

4 critical radii 

SCISi/C+O 

0.5 

20 

- 

- 

- 

Si/C+O interface 

SCIC+O/He 

0.5 

20 

- 

- 

- 

C+O/He interface 

SCIHe/H 

0.5 

20 

- 

- 

- 

He/H interface 

SCIC+O/He/H 

0.5 

20 

- 

- 

- 

C+O/He and He/H interfaces 

SC2p50m8f5 

0.5 

8 

5 

- 

- 

radial & angular modulated 

SC2p50m20f5 

0.5 

20 

5 

- 

- 

radial & angular modulated 

SC2pl0m8f5 

0.1 

8 

5 

- 

- 

radial & angular modulated 

SC2pl0m20f5 

0.1 

20 

5 

- 

- 

radial & angular modulated 

SC2p50m20fl 

0.5 

20 

1 

- 

- 

radial & angular modulated 

SC3p50fl 

0.5 

- 

1 

- 

- 

radial modulated 

SC3p50f5 

0.5 

- 

5 

- 

- 

radial modulated 

Bipo+Clp50m20 

0.5 

20 

- 

4.0 

- 

4 critical radii 

AsyEqua+Clp50m20 

0.5 

20 

- 

- 

1.8 

4 critical radii 

Asph+Clp50m20 

0.5 

20 

- 

4.0 

1.8 

4 critical radii 

Asph+Clp25m20 

0.25 

20 

- 

4.0 

1.8 

4 critical radii 

Asph+Clpl0m20 

0.10 

20 

- 

4.0 

1.8 

4 critical radii 

Asph+Clp50m20(*) 

0.5 

20 

- 

4.0 

1.8 

4 critical radii 

Asph+Clp25m20(**) 

0.25 

20 

- 

4.0 

1.8 

4 critical radii 


Note. — SCI, SC2, and SC3 are three spherical SN explosion cases presented in Section 2.3.1. Bipolar explosion is 
labeled as “Bipo” and equatorially asymmetric explosion is labeled as “AsyEqua”. The combination of bipolar and equatorially 
asymmetric explosions is labeled as “Asph”. The four critical radii in the CCSN progenitor are 3 x 10 8 cm at Si/C+O interface, 
3 X 10 9 cm at oxygen shell burning position, 6 x 10 9 cm at C+O/He interface, and 5 x 10 10 cm at He/H interface, respectively, 
eo is the density perturbation amplitude, m is the wave number of perturbations in the angular direction, and / is the scaling 
factor of the perturbation length scale. The ratio of the initial radial velocity along the polar axis to that along the equatorial 
axis Vp Q \/Ueq is given for bipolar explosion models. The ratio of the initial radial velocity v U p/v^ own at 6 = 0 to that at 6 = 7v 
is given for equatorially asymmetric models. All models have input energies of 2.5 X 10 51 erg; exceptions are those labeled by 
in which we set the input energy to be 1.7 X 10 51 erg and the one labeled by “**” in which we set the input energy to be 
2.0 x 10 51 erg. The model of SClp50m8(O) has a set of random numbers, which are different from the model of sclp50m8. 





Tabic 2: Main results of simulation models 


Model/Result 

Time 

Explosion 

“Ni tot 

56 Ni rati o 

44 Ti to t 

..max 

u r 

FigNo. 

unit 

(s) 

(10 51 erg) 

(M 0 ) 

(%) 

(io- 4 m 0 ) 

(km s" 1 ) 


Fiducialmodel 

6959 

1.91 

0.21 

- 

2.8 

1700 

- 

Basicmodel 

6982 

1.87 

0.19 

- 

3.5 

2000 

- 

SClp50m8 

6577 

1.80 

0.14 

23.1 

2.7 

4800 

Fig. 3 

SClp50m8(O) 

6557 

1.83 

0.15 

21.7 

2.1 

3700 

Fig. 4 

SClp50m20 

7062 

1.84 

0.17 

5.6 

7.2 

4200 

Fig. 5 

SClp50m20(*) 

9006 

1.08 

0.13 

0.9 

3.2 

4100 

Fig. 13 

SClpl0m8 

6978 

1.89 

0.19 

- 

2.1 

1700 

- 

SClpl0m20 

7029 

1.89 

0.20 

- 

2.8 

1500 

Fig. 6 

SCISi/C+O 

7005 

1.87 

0.18 

- 

5.5 

1700 

- 

SCIC+O/He 

7013 

1.86 

0.16 

0.7 

2.7 

3100 

- 

SCIHe/H 

7196 

1.88 

0.19 

0.6 

2.6 

2100 

- 

SCIC+O/He/H 

7087 

1.83 

0.15 

11.5 

2.3 

4800 

Fig. 7 

SC2p50m8f5 

5785 

2.01 

0.18 

8.3 

3.0 

3700 

Fig. 8 

SC2p50m20f5 

6618 

1.96 

0.21 

5.1 

3.8 

4200 

Fig. 9 

SC2pl0m8f5 

6762 

1.94 

0.21 

- 

3.0 

1800 

- 

SC2pl0m20f5 

6892 

1.93 

0.21 

- 

3.0 

1600 

Fig. 10 

SC2p50m20fl 

6636 

1.89 

0.21 

0.1 

4.4 

2800 

Fig. 11 

SC3p50fl 

6964 

1.88 

0.20 

- 

3.2 

1900 

Fig. 12 

SC3p50f5 

6330 

2.01 

0.19 

- 

2.9 

2000 

- 

Bipo+SClp50m20 

6714 

1.87 

0.09 

10.6 

6.9 

6000 

Fig. 14 

AsyEqua+SClp50m20 

6365 

1.88 

0.12 

7.0 

9.4 

6000 

Fig. 15 

Asph+SClp50m20 

6136 

1.85 

0.09 

19.4 

7.5 

5300 

- 

Asph+SClp25m20 

6101 

1.91 

0.09 

9.2 

4.0 

4200 

- 

Asph+SClpl0m20 

6382 

1.94 

0.09 

1.4 

2.7 

2100 

- 

Asph+SClp50m20(*) 

8133 

1.11 

0.06 

10.0 

4.9 

4500 

- 

Asph+SClp25m20(**) 

6984 

1.43 

0.07 

6.5 

2.8 

3400 

Fig. 16 


Note. — The physical time at which we stop the simulation is labeled as “Time” and the explosion energy 
is labeled as “Explosion”. The total masses of 56 Ni and 44 Ti are written as 56 Ni to t and 44 Ti to ti respectively. 
All of these element mass values are calculated from the 2D FLASH code in which the limited nucleosynthesis 
network is adopted. The ratio of the high velocity 56 Ni (> 3000 km s~ 4 ) mass to the total mass of 56 Ni 
is written as 56 Ni rat ; 0 . The maximum radial velocity of 56 Ni corresponding to the value of A M/M larger 
than 1.0 x 10 -3 in radial velocity profiles is labeled as u™ ax , where AM is the mass in the velocity range of 
v ~ v + Av, and M is the total mass. “FigNo.” is the figure number corresponding to the related simulation 
model in this paper. 



Table 3: Minimum and maximum radial velocities (u™ in and v™ ax ) of each element 
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